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Abstract 

We investigate geometric percolation and scaling relations in suspensions of nanorods, 
covering the entire range of aspect ratios from spheres to extremely slender nee¬ 
dles. A new version of connectedness percolation theory is introduced and tested 
against specialised Monte Carlo simulations. The theory accurately predicts perco¬ 
lation thresholds for aspect ratios of rod length to width as low as 10. The percola¬ 
tion threshold for rod-like particles of aspect ratios below 1000 deviates significantly 
from the inverse aspect ratio scaling prediction, thought to be valid in the limit of 
infinitely slender rods and often used as a rule of thumb for nano-fibres in composife 
maferials. Hence, mosf fibres fhaf are currenfly used as fillers in composife maferials 
cannof be regarded as practically infinilely slender for fhe purposes of percolafion 
fheory. Comparing percolafion fhresholds of hard rods and new benchmark resulfs 
for ideal rods, we find fhaf (i) for large aspecf rafios, fhey differ by a factor fhaf is 
inversely proporfional fo fhe connecfivify disfance befween fhe hard cores, and (ii) 
fhey approach fhe slender rod limif differenlly. 
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Connectivity percolation is the transition in which isolated clusters of solid par¬ 
ticles in a fluid (or of voids in a solid) become connected in some sense to form a 
system-spanning network. This network has a significant effect on the mechanical 
and transport properties of the material on a macroscopic scale iTj. If, for exam¬ 
ple, an electrically insulating polymer is mixed with conductive fibres such as car¬ 
bon nano tubes, the conductivity of the composite increases by ten or more orders of 
magnitude near the percolation transition of the filler material ll2j|. Given the tech¬ 
nological relevance to opto-electronics, photovoltaics and electromagnetic radiation 
shielding, it is no surprise that a large research effort is currently being invested to 
understand how the formulation and processing of a composite influence the perco¬ 
lation threshold as well as its physical properties beyond the threshold El]. 

The topic of percolation originates from studies on fluid flow in porous media, 
relevant, for example, to oil extraction. It has since been extensively studied theo¬ 
retically and computationally, both on and off lattice, the latter particularly (but not 
exclusively) for ideal, non-interacting bodies. An important scientific motivation for 
these studies is the critical behavior that the percolation transition shares with phase 
transitions [4]. For spherical particles, the impact of repulsive and attractive inter¬ 
actions on continuum percolation has received considerable attention [T|, while for 
non-spherical particles such as nano-wires, current understanding is much sketchier, 
despite their industrial interest as fillers in composite materials. 

Fibre-like fillers have been modelled as cylinders, spherocylinders and ellipsoids 
in theoretical studies B 0 H11 [la [m, in simulations where interactions are 
ignored ifT^ [131 [Ml [El, and in simulations where the particles interact via a hard 
excluded volume HU [m [El HI- However, there is no systematic test of theory 
against simulation over a large range of aspect ratios for interacting particles. In 
particular, the intermediate regime, in which the length is tens or a few hundred 
times the thickness, has not been addressed yet, even though most fibres used in 
realistic materials fall into this range. The reasons for this are that it is very time- 
consuming to simulate interacting particles of high length-to-width aspect ratios, and 
that analytical theories that are thought to be accurate in the limit of infinite aspect 
ratio are difficult to extend to finite values |[20ll . 

In this letter we present a combined theoretical and simulation study of perco¬ 
lation in suspensions of hard spherocylinders that spans an unprecedented range of 
aspect ratios from spheres to slender rods of aspect ratio up to 1000. We show that 
hard-core interactions (i) change the approach to the theoretically expected scaling of 
the percolation threshold with the inverse aspect ratio, (ii) shift the percolation thresh¬ 
old to larger values than that for ideal rods by a factor that converges only slowly with 
respect to rod length and (iii) cause the core packing fraction at percolation to exhibit 
a maximum for small aspect ratios. These results show that the ideal (penetrable) par¬ 
ticle model has limited predictive value for actual rod systems and that even for very 
long rods, a finite-length correction is necessary to get quantitative predictions. We 
show that this correction can be obtained explicitly from connectedness percolation 
theory, using a sensible ansatz for the connectedness direct correlation function. 

Before discussing the results, we briefly describe our simulation method and 
theory. We have generated configurations of hard spherocylinders at fixed particle 
number N and volume V, using cubic simulation cells of length The 

spherocylinders consist of a cylinder of length L and diameter D, capped with hemi¬ 
spheres of the same diameter. Hence, the surface of the spherocylinder consists of 
all points lying a distance D/2 from a line segment of length L. The full aspect ratio 
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Figure 1: Schematic illustration of two spherocylinders with overlapping contact shells. L 
is the length of the line segment (white) of each spherocylinder, which is common to the 
hard core and to the contact shell, D is the diameter of the hard core (dark), and X is the 
diameter of the contact shell (light). 


of a hard spherocylinder, including the caps, is L/D + 1. The core is strictly im¬ 
penetrable but is surrounded by a notional contact shell that is used to define when 
two spherocylinders are considered to be connected. The surface of the contact shell 
is a spherocylinder that shares the same line segment as the core but has diameter 
X instead of D. Hence, the surface of the contact shell lies at a uniform distance 
{X — D)/2 from the surface of the core (see Fig.[T]l. The full aspect ratio of the con¬ 
tact shell is therefore L/X+l. Two hard spherocylinders are connected if their shells 
overlap, and clusters are defined by contiguous pairwise connections. 

We sample the fraction P{(^-,Lx) of independent configurations that contain a 
percolating cluster as a function of the packing fraction (j) of the cores for a given 
cell length Lx. The packing fraction is defined by (|) = Nvcore/V, where Vcore = 
hD^ (f+ f) is the volume of one spherocylinder. To detect a percolating cluster 
in a cubic simulation cell with periodic boundary conditions, we require that a clus¬ 
ter must connect periodic images of its constituent particles in at least one of the 
periodic directions. This “wrapping” criterion is somewhat more costly to evaluate 
than the simpler “spanning” criterion, which only requires that a cluster connects two 
opposite boundaries of the simulation cell. However, wrapping clusters are a more 
accurate representation of a percolating cluster in the macroscopic limit because such 
clusters are infinite when the simulation cell is replicated through its boundary con¬ 
ditions. In contrast, spanning clusters merely form an array of large but disconnected 
clusters when the cell is replicated. Furthermore, wrapping probabilities follow uni¬ 
versal scaling functions |[^ —a feature that we will exploit in the simulations of 
ideal rods to mitigate the small systematic errors arising from the inevitably finite 
size of the simulation cell. At finite Lx, P{(^\Lx) is a sigmoidal function of (|), be¬ 
coming a sharp step function as Lx —)• 0 °. The curves for different Lx have a common 
crossing point ll^ . typically just below P = 0.5. To make the simulations tractable, 
we identify the hard-rod percolation threshold (|)p as the point /’((|)p;Ljc) = 0.5 using 
a cell length up to Lx = lOL (never less than Lx = 15D). Fig. |2] shows examples of 
percolating clusters for L/D = 200 and L/D = 10. 

Properly equilibrated configurations of impenetrable rods cannot be generated 
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Figure 2: Simulation snapshots of percolating clusters. The rods are colour coded accord¬ 
ing to their orientation. Left panel: L/D = 200, right panel: L/D = 10. 

by sequential random insertion If22l . Configurations must therefore be obtained by 
Monte Carlo displacements and rotations of particles, rejecting any trial move that 
would generate an overlap of cores and accepting all others. This procedure becomes 
computationally costly for slender rods, where the simulation cell and number of par¬ 
ticles must be large to allow for an accurate computation of the percolation threshold. 

To find overlaps quickly, we use a mefhod in which fhe spherocylinders are nofion- 
ally divided info small segmenfs so fhaf the overlap detection need only be performed 
for neighboring segments i2^ . This method scales linearly with particle number but 
requires a large amount of memory. 

We have also calculated the percolation threshold for fully penetrable (ideal) 
spherocylinders, consisting of the contact shell with no hard core. By definition, 

(j) = 0 for ideal rods, but hard and ideal spherocylinders may be compared using 
a packing fraction renormalised to the volume of the contact shell, r] = A^Vsheii/T, 
where Vsheii = + the hypothetical volume fraction of the shells in the 

absence of shell overlaps. Independent equilibrium configurations of ideal sphero- 
cylinders are easily generated by placing rods with random positions and orientations 
within the simulation cell. For ideal rods we used cells lengths up to = 6L (never 
less than 16X) and identihed the percolation threshold from the crossing point of 
P{r\',Lx) at two cell sizes L^, sampling 20000 conhgurations at each r] and L^. This 
procedure minimises the systematic errors due to finite cell sizes and reduces statis¬ 
tical uncertainty in the ideal-rod percolation thresholds to 0.1%. 

Our theoretical predictions are based on connectedness percolation theory. The 
percolation threshold is defined as fhe filler frachon (|) af which fhe mean clusfer size 
diverges and is equal fo Vcore/ {{C^))' 123 ■ Here, Vcore is again the particle volume, 

C+ the spatial Fourier transform of the connectedness direct correlation function, 
C'^(r,u,u'), at zero wave vector. The direct correlation function is a function of 
the vector r connecting the centres of mass of the particles, and their main-body- 
axis vectors u and u'. The angular brackets denote an orientational averaging over 
these two directors. In the limit of long thin rods the second virial approximation 
is accurate ll20l . Within the second virial approximation, C+(r,u,u') = /+(r,u,u'), 
where /+ = exp(— 17+) is the connectedness Mayer function with 17+ = 17+(r,u,u') 
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the connectedness interaction potential scaled to the thermal energy. 17+ = oo for 
non-connected configurations and 17+ = 0 for connecting ones, i.e., configurations 
for which the connectivity shells overlap but the hard cores do not. 

To go beyond the second virial approximation, we invoke a Lee-Parsons type of 
approximation that has proven remarkably accurate in predicting the phase behav¬ 
ior of hard rods and mixtures of hard rods and hard spheres flM . It is based on an 
interpolation between the Percus-Yevick equation of state for hard spheres and the 
second virial equation of state for hard rods. In the context of connectedness perco¬ 
lation, it can be written as C+ = /+ x (1 — 3(|)/4)/( 1 — (|))^, given the known relations 
between the direct correlation function and the connectedness variant of it |[20ll25l . 
This then gives for the percolation threshold an explicit expression in terms of the 
ratio Y= L/D and the dimensionless measure of the connectivity range a = X/D — 1, 


2(l+^)-2(l + l^)^ 
3(l-hf^) 


( 1 ) 


where 


^(Y,a) = 


1 + ^^U 

^ + 3y; y 
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l-foc) ~l)+^((l + ct) —l)-|-OC 


( 2 ) 


To obtain numerical results for the percolation threshold, we must select a value 
for the connectivity criterion X. The separation X — D of the core spherocylinder sur¬ 
faces at the connectivity cut-off (see Fig. [U should be characteristic of the distance 
over which electron tunnelling between nanorods decays ll26l . This distance depends 
on the details of both the nanorods and the medium in which they are suspended Q 
and a full quantum mechanical treatment is a formidable task. However, if the en¬ 
ergetic barrier AE for tunnelling can be measured or calculated then an estimate for 
X — D can be obtained from the tunnelling length (/i^/2/7ieA£')^/^ through a rectan¬ 
gular barrier, where nig is the electron mass. For suspensions of carbon nanotubes 
(diameters on the order of a nanometre), the tunnelling length typically lies at a frac¬ 
tion of a nanometre [3]]. We have therefore taken X/D = 1.2 as a representative value 
for most of our results, but will consider the effects of altering this value later. 

In Fig. [3] we present our simulation results (circles) for the percolation threshold 
of hard rods in terms of the renormalised volume fraction r] versus the full aspect 
ratio L/X+l, from spheres (aspect ratio 1, L/D = 0) to very slender rods with core 
L/D = 1000. For comparison, we also show (i) our simulation results for ideal rods 
(pluses), (ii) the phenomenological expression for ideal cylinders of the same aspect 
ratio, obtained by Mutiso et al. using a fit to simulation data on somewhat shorter rods 
HI (dashed line), and (iii) our theoretical prediction (solid line). Our theory agrees 
quantitatively with the simulations of hard rods for aspect ratios above 10, and semi- 
quantitatively below that. The largest discrepancy between simulation and theory 
occurs for aspect ratios approaching zero and amounts to 25 % for our choice of 
connectivity range. This is to be expected given the level of accuracy of the Percus- 
Yevick approximation for the percolation of hard spheres lITTl . The Mutiso fit for 
ideal cylinders ifTSTl captures the qualitative deviations from the long-rod limit, but 
only agrees with our simulations within statistical error over a narrow middle range 
of aspect ratios. 

A number of conclusions can be drawn from Fig. [3] First, we may compare the 
accurate results of the Monte Carlo simulations with the asymptotic scalings for the 
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Figure 3: Percolation threshold Tjp as a function of aspect ratio for ideal rods and hard 
rods. The approach to the scaling regime, where (L/X+ l)ri would be constant in each 
case, is highlighted in the inset. Statistical errors in the ideal-rod results are smaller than 
the symbols. Tabulated data for the main plot are available as supplementary material. 
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Figure 4: Percolation threshold as a function of the contact shell thickness, i.e., the 
surface-to-surface distance criterion for connectivity. 

percolation thresholds of long rods. The long-rod scaling predictions are rjp ~ \j2L 
for ideal rods 12^ and rjp ~ jlLQ^ — D') for hard rods. Both these results can be 
obtained from connectedness percolation theory within the second virial approxima¬ 
tion in the limit L D Q. However, the same prediction for ideal rods was also 
made in earlier work 1^12^ . using a conjecture based directly on average excluded 
volume and the number of contacts between objects at the percolation threshold. In 
both the ideal- and hard-rod cases, the product (L/k-|- Ijr] should approach a con¬ 
stant value with increasing rod length. In the inset of Fig. [3l however, we see that 
the asymptotic scalings for the percolation threshold of long rods are only reached 
for aspect ratios in excess of several hundred 1201 . For hard rods, a constant value of 
(L/X -I- l)ri is reached slowly from below, while for ideal rods the respective plateau 
is approached from above after initially overshooting. 

Second, hard-core interactions seem to have a larger impact on the percolation 
threshold for large aspect ratios than for smaller ones. This observation agrees with 
previous Monte Carlo simulations of hard spheres ll^l22l . The smaller the connec¬ 
tivity range, the larger the difference between the percolation thresholds of hard and 
ideal rods. Fig. |4] illustrates this point, showing a significant impact of the connec¬ 
tivity range on the percolation threshold even for relatively short rods of aspect ratio 
10 and 20, itself quite accurately predicted by our theory. In the slender rod limit, 
the ratio of percolation thresholds of ideal to hard rods should be proportional to 
the thickness of the connectivity shell around the core, as (X —D)/k. The results in 
Fig. |3]bear out this theoretical prediction but, again, the limiting behavior is reached 
for the very longest of the rods included in our simulations. 

In Fig. |5]we focus on very short hard rods (1.3 < L/D -\-1 < 3). The connectiv¬ 
ity criterion was again set to k/D = 1.2. For short rods, if we define the percolation 
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Figure 5: Percolation threshold for short hard rods with connectivity criterion X/D = 1.2 
from simulation. If plotted as a function of the physical volume fraction of the hard cores 
(left-hand panel), the percolation threshold is non-monotonic in the aspect ratio, while the 
notional volume fraction of the connectivity shells (right-hand panel) shows no extremum. 


threshold in terms of the physical volume fraction (|)p due to the cores of the particles, 
we find the threshold to be a non-monotonic function of the aspect ratio (left-hand 
panel of Fig. [5]), with a maximum close to L/D = 0.7 (full core aspect ratio 1.7). Our 
theory also predicts a maximum at short aspect ratios. Maximisation of (|)p in Eq. O 
with respect toy = L/D gives L/D = |(—1-|- yjl + \l(X + 6o?) for the aspect ratio 
at maximum percolation threshold. Using the connectivity range a = X/D — 1 = 0.2 
from the simulations, the theoretical expression evaluates to L/D k, 1.4, i.e., about 
a factor ot 2 larger than the simulation result. One might be tempted to speculate 
that the maximum in (j)p is related to the fact that slightly aspherical objects, such as 
ellipsoids, pack more efficiently than spheres |[^ . However, if the data are plotted 
in terms of the renormalised volume fraction rjp, the curve becomes monotonic (see 
Fig. [21 right-hand panel). Thus the maximum arises from the definition of connectiv¬ 
ity. Crucially, however, it is (|) that corresponds to the experimental volume fraction 
of rod-like filler. Hence, real composite materials should indeed show a local maxi¬ 
mum in loading a the percolation threshold when short filler particles are used. The 
renormalised volume fraction r] that includes the effective, penetrable contact shell 
is generally not accurately known a priori. 

In summary, we have introduced a new version of connectedness percolation the¬ 
ory which accurately predicts percolation thresholds for nanorods over a large range 
of aspect ratios. We show that for aspect ratios below 1000 the percolation threshold 
deviates significantly from the inverse aspect ratio scaling prediction. Hence, caution 
is required when making predictions about systems in the intermediate regime of as¬ 
pect ratios that is relevant to typical materials applications. We have also presented 
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new simulation data for hard and ideal rods, whieh are the first to eover aspeet ratios 
from a sphere to very slender rods, providing a new benehmark. Ideal rods—a eonve- 
niently simple model—differ both quantitatively and qualitatively from more phys- 
ieal models that inelude an impenetrable eore, and the model’s limitations should 
therefore not be ignored in the eontext of real nanorod systems. Hard-eore inter- 
aetions ehange the approaeh to the theoretieally expeeted sealing of the pereolation 
threshold with the inverse aspeet ratio and shift the pereolation threshold to larger 
values than that for ideal rods by a faetor that eonverges only slowly with respeet to 
rod length. We have shown that the eorreetion to the long-rod sealing regime that is 
needed at all realistie aspeet ratios of hard nanorods ean be obtained from eonneet- 
edness pereolation theory, demonstrating the strength of this approaeh. 
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Supplemental Table 1: Percolation threshold of hard spherocylinders with cylin¬ 
drical portion of length L, diameter D and contact shell thickness X = 1 .2D, from 
Monte Carlo simulations (blue circles in Figure 2 of the article). 


core aspect ratio 
L/D 

overall aspect ratio 
L/X+l 

percolation threshold 
Tip 

uncertainty 

5rip 

0.3 

1.25 

0.349 

0.002 

0.5 

1.416667 

0.347 

0.002 

0.6 

1.5 

0.346 

0.002 

0.7 

1.583333 

0.343 

0.002 

0.732051 

1.610042 

0.343 

0.002 

0.8 

1.666667 

0.342 

0.002 

1.0 

1.833333 

0.337 

0.002 

2.0 

2.666667 

0.311 

0.002 

10.0 

9.333333 

0.17974 

0.00009 

20.0 

17.666667 

0.11848 

0.00003 

30.0 

26 

0.08892 

0.00002 

50.0 

42.666667 

0.05930 

0.00007 

70.0 

59.333333 

0.0447 

0.0003 

100.0 

84.333333 

0.03250 

0.00004 

200.0 

167.666667 

0.0171 

0.0004 

1000.0 

834.333333 

0.00359 

0.00006 
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Supplemental Table 2: Percolation threshold of ideal spherocylinders with cylin¬ 
drical portion of length L and diameter X from Monte Carlo simulations (red pluses 
in Figure 2 of the article). The fractional uncertainty in the percolation threshold is 
5t1p/t1p = 0.001 


overall aspect ratio 
L/X+l 

percolation threshold 
Tip 

1 

0.3424 

1.1 

0.3403 

1.2 

0.3373 

1.3 

0.3320 

1.4 

0.3256 

1.5 

0.3188 

1.6 

0.3122 

1.8 

0.2969 

2 

0.2823 

2.2 

0.2687 

2.5 

0.2495 

3 

0.2215 

4 

0.1772 

5 

0.1464 

7 

0.1072 

9 

0.08368 

11 

0.06814 

16 

0.04564 

21 

0.03389 

31 

0.02208 

41 

0.01619 

51 

0.01279 

61 

0.01044 

81 

0.007636 

101 

0.005992 

151 

0.003878 

201 

0.002846 

401 

0.001366 

601 

0.0008955 

1001 

0.0005281 

2001 

0.0002590 
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